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We modify the kinetic Ising model with Metropolis dynamics, allowing each spin to interact only 
with q spins randomly chosen from the whole system, which corresponds to the topology of a com¬ 
plete graph. We show that the model with q > 3 exhibits a phase transition between ferromagnetic 
and paramagnetic phases at temperature T*, which linearly increases with q. Moreover, we show 
that for q = 3 the phase transition is continuous and discontinuous for larger values of q. For q > 3 
the hysteresis exhibits oscillatory behavior - expanding for even values of q and shrinking for odd 
values of q. If only simulation results were taken into account, this phenomenon could be mistakenly 
interpreted as switching from discontinuous to continuous phase transitions or even as evidence of 
the so-called mixed phase transitions. Due to the mean-field like nature of the model we are able 
to calculate analytically not only the stationary value of the order parameter but also precisely 
determine the hysteresis and the effective potential showing stable, unstable and metastable steady 
states. The main message is that in case of non-equilibrium systems the hysteresis can behave in 
an odd way and computer simulations alone may mistakenly lead to incorrect conclusions. 


I. INTRODUCTION 


II. THE MODEL 


According to the modern theory of phase transitions, 
each phase transition can be described by an order pa¬ 
rameter, having a non-zero value in the ordered phase 
and vanishing in the disordered phase dHi. The be¬ 
havior of an order parameter allows to classify each 
phase transition as a continuous or discontinuous (first 
order). However, it has been noticed in a number of 
cases that the dichotomy between continuous and discon¬ 
tinuous transitions fails, in the sense that a jump of the 
order parameter coincides with power-law singularities 
[3-[l| or even with the absence of fundamental indicators 
of the first-order phase transitions, such as the hysteresis, 
metastable states and phase coexistence It can also 
happen that the phase transition is weakly discontinu¬ 
ous, i.e. the jump of the order parameter is small and, 
therefore, to decide on the type of the phase transition 
is quite difficult in computer simulations. In such a situ¬ 
ation measuring the hysteresis of the order parameter is 
a demanding task Q. 


In this paper we show that even measuring the hys¬ 
teresis based solely on computer simulations can be mis¬ 
leading. We propose here to modify the kinetic Ising 
model with Metropolis dynamics @,0 on a complete 
graph by assuming that in each elementary time step a 
randomly chosen spin interacts only with its q neighbors. 
We show that in such a model the hysteresis exhibits 
oscillatory behavior, expanding for even values of q and 
shrinking for odd values of q, which could be mistakenly 
interpreted as switching from discontinuous to continu¬ 
ous phase transitions or even as evidence of the so-called 
mixed phase transitions Hi &[1|- However, due to the 
mean-field like nature of the model we are able to calcu¬ 
late analytically not only the stationary value of an order 
parameter but also to determine the hysteresis and the ef¬ 
fective potential showing stable, unstable and metastable 
steady states, analogously as it was done for the q-voter 
model with noise 111. 


The idea to consider exactly q neighbors, no matter 
what is the actual number of neighbors on a given graph, 
is borrowed from the q-voter model originally pro¬ 

posed to introduce non-linearity in the voter dynamics 
at the microscopic level. Within the q-voter model, each 
spin is described by a dynamical variable S'i = ±1 and in¬ 
teracts with a set of q neighbors. If all q neighbors share 
the same state, the spin conforms to this state. In the 
other case the spin flips with probability e. It is worth 
to notice here, that the one-dimensional q-voter model 
with q = 2 is identical to the Ising model with general¬ 
ized zero-temperature Glauber dynamics in which a 
spin flips with probability p = 1 in the case of energy de¬ 
crease and with probability p = Wq in the case of energy 
conservation. If we denote Wq = e, the time evolution of 
a single spin for both models can be written as follows: 

( 1 with p = 1 if Si-i = 5*+! = 1, 

S[ = < -S^ with p = Wo if S'j-iS'j+i = -1, (1) 

y -1 with p = 1 if S^-i = Sj+i = -1, 

where for brevity we use the notation S'' = Si{t+At) and 

Si = Si{t) and Wq = 1 corresponds to the Metropolis, 
whereas Wq = 1/2 to the original Glauber dynamics 0. 

For higher dimensions, both models are not equivalent 
even in zero temperature and even for q equal to the num¬ 
ber of the nearest neighbors. The q-voter model requires 
a unanimous state of all q neighbors to influence spin Si, 
whereas for the Ising model a majority is sufficient, which 
follows from the Hamiltonian: 

H = -J2s.S,. ( 2 ) 

i,3 

However, one could consider the q-voter model with 
threshold r = 1/2 (i.e. majority needed to influence the 
spin) [T^ and then again both models would be equiva¬ 
lent at zero temperature. 
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The behavior of the Ising model described by the 
Hamiltonian in Eq. ([2]) under zero-temparature Glauber 
dynamics is very interesting, exhibiting a slow relax¬ 
ation related to a metastable state d, [ij, [T^. However, 
here we focus on another problem related to the kinetic 
Ising model, inspired by the analogy between the one¬ 
dimensional kinetic Ising model with zero-temperature 
Glauber dynamics and the g-voter model with q = 2. 
We ask the following question: what would be the behav¬ 
ior of a modified kinetic Ising model - in which every 
spin interacts with a set of q neighbors randomly eho- 
sen from the set of all its neighbors - if we introduced a 
temperature-like parameter T > 0? 

The algorithm of a single time step of the g-neighbor 
Ising model consists of 3 consecutive steps: 

1. Randomly choose a spin, Si, and from all its neigh¬ 
bors choose a subset of g neighbors, nuq. 

2. Calculate the value of the following function, based 
on the Hamiltonian in Eq. for the original state 
of the i-th spin: 

EiS.) = -S. Y. ( 3 ) 

jGrin, 

and the value of the same function for the flipped 
j-th spin, i.e. E(—Si). 

3. Flip the i-th spin with probability min[l, 
where AE = E{—Si) — E{Si). 

We would like to stress here, that minimizing the func¬ 
tion given in Eq. ([3]) does not necessarily lead to the 
minimization of the whole energy of the system given by 
Eq. ([2]), in contrast to the equilibrium Ising model in 
which we sum interactions over all nearest neighbors. As 
usual we choose the magnetization: 

1 ^ 

w(t) = (4) 

' i=l 

as an order parameter, which in the case of a complete 
graph fully describes the state of the system. 

III. RESULTS 

We investigate the model using an analytical approach 
and Monte Carlo simulations. The latter start from two 
types of initial conditions - fully ordered (m = 1), which 
corresponds to zero temperature, and completely random 
(to = 0), which corresponds to a high temperature. For 
each value of the temperature-like parameter T, we mea¬ 
sure the stationary value of the magnetization defined 
byEq. 0. We have checked that averaging over time 
gives the same result as averaging over samples. However, 
when using the time average it is easier to distinguish 
between continuous and discontinuous phase transitions 


q=3 q=4 




FIG. 1. Dependencies between steady values of the magneti¬ 
zation m and the temperature-like parameter T for 4 values 
of q. Symbols represent results obtained from Monte Carlo 
simulations for a system of size N = 10® and two types of 
initial conditions -- fully ordered (o) and disordered (*). Nu¬ 
merical results, see Eq. 0, are presented by lines - solid for 
stable solutions and dotted for unstable. 


looking solely at the order parameter as a function of T. 
Starting from two types of initial conditions should allow 
us to identify the type of the transition on the basis of 
the hysteresis, but as we will see later, the hysteresis can 
be also misleading. 

Dependencies between steady values of the magneti¬ 
zation TO and the temperature-like parameter T are pre¬ 
sented in Fig. [TJ A phase transition between ordered 
and disordered phases is observed for all values of g > 3. 
For g = 3 there is no jump in the order parameter and no 
hysteresis, which indicates a continuous phase transition. 
For g = 4 and g = 6 the jump of the order parameter and 
the hysteresis indicate a discontinuous phase transition. 
However, in the case of g = 5, we are not able to distin¬ 
guish between the continuous and discontinuous phase 
transition - the jump of the order parameter is observed 
and simultaneously there is no hysteresis, similarly as 
observed in 

Fortunately, thanks to the mean-field type of the model 
we are able, following the reasoning presented in PH, 
to write down equations that allow to calculate the sta¬ 
tionary value of the order parameter, as well as an effec¬ 
tive potential. The latter helps to distinguish between 
continuous and discontinuous phase transitions. 

Although the magnetization to is an order parameter 
of the system and all results will be presented in terms 
of TO, calculations are easier if we use the concentration 
of ‘up-spins’ which is equivalent to the probability that 
a randomly chosen spin is ‘up’. In a single time step 
At, three events are possible - the concentration of ‘up- 
spins’ c increases by 1/A, decreases by 1/A or remains 
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constant: 


7 +(c, T, q) = Prob {c —>■ c + An} , 

'y~{c,T,q) = Prob{c^ c - An] , 

l°{c,T,q) = Prob{c c} = 1 -7'^(c) - 7”(c), (5) 


where: 


Q 


7'lc,T,q) = ^(^)c«-'=(l-c)'=+i 


7 (c,T,g) = ^ L )(l-c) 


A:=0 


min 


1 , 

I g#(9-2fc) 


( 6 ) 


In the stationary state we expect that the probability 
of growth 7 + (c, T, q) is equal to the probability of loss 
'y~{c,T,q). Therefore: 


F{c, T, q) = 7 + (c, T,q)-^- (c, T,q)=0, (7) 


where F{c, T, q) can be treated as an effective force that 
drives the concentration c up or down m- Solving an¬ 
alytically Eq. 0, i.e. finding the stationary value of c 
as a function of T for an arbitrary value of q, is impossi¬ 
ble, but we can easily do it numerically. The results are 
denoted by lines in Fig. [TJ 

Having the effective force F{c,T,q) we can also calcu¬ 
late the effective potential: 

Vic,T,q) = -J F{c,T,q)dc, (8) 

which, as seen in Fig. [2l allows to distinguish between 
stable (minima of the potential) and unstable (maxima 
of the potential) solutions of Eq. 0. Furthermore, it 
allows to distinguish between continuous and discontin¬ 
uous phase transitions, in the latter case showing phase 
coexistence and metastable states. Finally, it allows to 
determine the transition point which in case of a discon¬ 
tinuous phase transition coincides with the value of T, for 
which minima corresponding to disordered and ordered 
phases are equal. 

Having the effective force 0 and the effective potential 
([8|) we can calculate numerically not only the stationary 
value of the order parameter 'rn{T) for an arbitrary value 
of q, but also determine the transition temperature T* 
and the width of the hysteresis defined here as the dis¬ 
tance between the spinodal lines (see Fig. 0. 

More precisely, the latter can be calculated from the 
potential 0. For low values of T there are two minima 
that correspond to ordered phases. Then at T = Ti the 
third minimum appears; it corresponds to the disordered 
phase but it shallower than other two, i.e. the disordered 
state is metastable. At T = T* all three minima are 
equal, which corresponds to the transition point. Above 
this value the middle minimum corresponding to the dis¬ 
ordered phase is the deepest and the other two represent 
metastable ordered states. Finally, above T = T 2 there 
is only one minimum - the disordered state is the only 



FIG. 2. Potentials given by Eq. m for q — S (top panels) 
and q = b (bottom panels). Left panels represent potentials 
just below and right panels just above the transition point. 
It is seen that for q = S there are two ordered phases below 
the transition point and a single disordered phase above the 
transition point, which is typical for continuous phase transi¬ 
tions. For q = 5 the phase coexistence and metastable states 
are seen, which indicates a discontinuous phase transition. 
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FIG. 3. Dependence between the phase transition value of T 
(i.e. r*; left panel) or the width of the hysteresis (i.e. AT; 
right panel) and the number of neighbors (i.e. q). 


possible state of the system. The distance between the 
spinodal lines AT = T 2 — T 1 determines the width of the 
hysteresis and is presented in the right panel of Fig. [3] 

In Fig. [3] it can be seen that the transition tempera¬ 
ture T* increases linearly with q, whereas the hysteresis 
exhibits an oscillatory behavior, expanding for even val¬ 
ues of q and shrinking for odd values of q. Interestingly, 
for g = 5 the width of the hysteresis is roughly zero and 
even for g = 7 it is difficult to see any hysteresis in com¬ 
puter simulations. This fact has initially led us to the 
wrong conclusion of oscillatory switching from continu¬ 
ous to discontinuous phase transitions. 























4 


IV. CONCLUSIONS 


It is believed that the hysteresis is the main indicator 
of discontinuous phase transitions ii- However, as we 
have shown here, the behavior of the hysteresis might 
be quite unexpected. In the simple model studied here, 
we have observed a jump of the order parameter and si¬ 
multaneously no hysteresis for q = 5. Even for larger 
odd values oi q {q = 7,9) the hysteresis was invisible in 
computer simulations. If the model was not solvable ana¬ 
lytically, the results of computer simulations could drive 
us to the wrong conclusion of a hybrid phase transition in 
which the jump of the ordered parameter coincides with 
no hysteresis. Luckily, due to the mean-field character of 
the model, we could distinguish between the two types 
of transitions on the basis of the effective potential. 

In equilibrium statistical mechanics, it is common that 
systems, which exhibit a discontinuous phase transition 
in high space dimensions, may display a continuous tran¬ 


sition below a certain upper critical dimension Q. We 
expect that at the same time the hysteresis monotoni- 
cally decays, reaching zero at the upper critical dimen¬ 
sion. Here we have similar situation, despite the fact that 
we do not change the dimension but only the number of 
neighbors q. For q > 3 the transition is discontinuous, 
for g = 3 continuous and for g < 3 there is no phase 
transition in the system. However, at the same time the 
hysteresis does not decay monotonically but oscillates. 
Therefore, the main message of our paper is a warning 
that in case of non-equilibrium systems the hysteresis can 
behave in an odd way and computer simulations alone 
may mistakenly lead to incorrect conclusions. 
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